home *** CD-ROM | disk | FTP | other *** search
/ SGI Developer Toolbox 6.1 / SGI Developer Toolbox 6.1 - Disc 4.iso / src / demos / VGX / shadows / trackball.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-08-02  |  5.2 KB  |  216 lines
  1. /*
  2.  * Copyright 1991, 1992, 1993, 1994, Silicon Graphics, Inc.
  3.  * All Rights Reserved.
  4.  *
  5.  * This is UNPUBLISHED PROPRIETARY SOURCE CODE of Silicon Graphics, Inc.;
  6.  * the contents of this file may not be disclosed to third parties, copied or
  7.  * duplicated in any form, in whole or in part, without the prior written
  8.  * permission of Silicon Graphics, Inc.
  9.  *
  10.  * RESTRICTED RIGHTS LEGEND:
  11.  * Use, duplication or disclosure by the Government is subject to restrictions
  12.  * as set forth in subdivision (c)(1)(ii) of the Rights in Technical Data
  13.  * and Computer Software clause at DFARS 252.227-7013, and/or in similar or
  14.  * successor clauses in the FAR, DOD or NASA FAR Supplement. Unpublished -
  15.  * rights reserved under the Copyright Laws of the United States.
  16.  */
  17. /*
  18.  *    Implementation of a virtual trackball.  See trackball.h for the
  19.  * interface to these routines.
  20.  *    Implemented by Gavin Bell, lots of ideas from Thant Tessman and
  21.  * the August '88 issue of Siggraph's "Computer Graphics," pp. 121-129.
  22.  */
  23. #include <stdio.h>
  24. #include <math.h>
  25. #include "trackball.h"
  26.  
  27. /* Size of virtual trackball, in relation to window size */
  28. float trackballsize = 0.4;
  29.  
  30. /* Function prototypes for local functions */
  31. float tb_project_to_sphere(float, float, float);
  32. void normalize_euler(float *);
  33.  
  34. /*
  35.  * Ok, simulate a track-ball.  Project the points onto the virtual
  36.  * trackball, then figure out the axis of rotation, which is the cross
  37.  * product of P1 P2 and O P1 (O is the center of the ball, 0,0,0)
  38.  * Note:  This is a deformed trackball-- is a trackball in the center,
  39.  * but is deformed into a hyperbolic solid of rotation away from the
  40.  * center.
  41.  * 
  42.  * It is assumed that the arguments to this routine are in the range
  43.  * (-1.0 ... 1.0)
  44.  */
  45. void
  46. trackball(float *e, float p1x, float p1y, float p2x, float p2y)
  47. {
  48.     float p1[3];
  49.     float p2[3];
  50.     float d[3];
  51.     float phi;    /* how much to rotate about axis */
  52.     float a[3];    /* Axis of rotation */
  53.  
  54.     vzero(a);
  55.  
  56.     if (p1x == p2x && p1y == p2y)
  57.     {
  58.         vzero(e);    /* Zero rotation */
  59.         e[3] = 1.0;
  60.         return ;
  61.     }
  62.  
  63. /*
  64.  * First, figure out z-coordinates for projection of P1 and P2 to
  65.  * deformed sphere
  66.  */
  67.     vset(p1, p1x, p1y, tb_project_to_sphere(trackballsize, p1x, p1y));
  68.     vset(p2, p2x, p2y, tb_project_to_sphere(trackballsize, p2x, p2y));
  69.  
  70. /*
  71.  *    Now, we want the cross product of P1 and P2
  72.  *  (Or cross product of p2 and p1... this was determined by trial and
  73.  * error, so there may be a compensating mathematical boo-boo
  74.  * somewhere else).
  75.  */
  76.     vcross(p2, p1, a);
  77.  
  78. /*
  79.  *    Figure out how much to rotate around that axis.
  80.  */
  81.     vsub(p1, p2, d);
  82.     phi = fsin(vlength(d) / (2.0 * trackballsize));
  83.  
  84.     axis_to_euler(a, phi, e);
  85. }
  86.  
  87. /*
  88.  *    Given an axis and angle, compute euler paramaters
  89.  */
  90. void
  91. axis_to_euler(float *a, float phi, float *e)
  92. {
  93.     vnormal(a);    /* Normalize axis */
  94.     vcopy(a, e);
  95.     vscale(e, fsin(phi/2.0));
  96.     e[3] = fcos(phi/2.0);
  97. }
  98.  
  99. /*
  100.  * Project an x,y pair onto a sphere of radius r OR a hyperbolic sheet
  101.  * if we are away from the center of the sphere.
  102.  */
  103. float
  104. tb_project_to_sphere(float r, float x, float y)
  105. {
  106.     float d, t, z;
  107.  
  108.     d = sqrt(x*x + y*y);
  109.     if (d < r*M_SQRT1_2)    /* Inside sphere */
  110.     {
  111.         z = sqrt(r*r - d*d);
  112.     }
  113.     else    /* On hyperbola */
  114.     {
  115.         t = r / M_SQRT2;
  116.         z = t*t / d;
  117.     }
  118.     return z;
  119. }
  120.  
  121. /*
  122.  *    Given two rotations, e1 and e2, expressed as Euler paramaters,
  123.  * figure out the equivalent single rotation and stuff it into dest.
  124.  * 
  125.  * This routine also normalizes the result every COUNT times it is
  126.  * called, to keep error from creeping in.
  127.  */
  128. #define COUNT 100
  129. void
  130. add_eulers(float *e1, float *e2, float *dest)
  131. {
  132.     static int count=0;
  133.     register int i;
  134.     float t1[3], t2[3], t3[3];
  135.     float tf[4];
  136.  
  137.     vcopy(e1, t1); vscale(t1, e2[3]);
  138.     vcopy(e2, t2); vscale(t2, e1[3]);
  139.     vcross(e2, e1, t3);
  140.     vadd(t1, t2, tf);
  141.     vadd(t3, tf, tf);
  142.     tf[3] = e1[3] * e2[3] - vdot(e1, e2);
  143.  
  144.     for (i = 0 ; i < 4 ;i++)
  145.     {
  146.         dest[i] = tf[i];
  147.     }
  148.  
  149.     if (++count > COUNT)
  150.     {
  151.         count = 0;
  152.         normalize_euler(dest);
  153.     }
  154. }
  155.  
  156. /*
  157.  * Euler paramaters always obey:  a^2 + b^2 + c^2 + d^2 = 1.0
  158.  * We'll normalize based on this formula.  Also, normalize greatest
  159.  * component, to avoid problems that occur when the component we're
  160.  * normalizing gets close to zero (and the other components may add up
  161.  * to more than 1.0 because of rounding error).
  162.  */
  163. void
  164. normalize_euler(float *e)
  165. {    /* Normalize result */
  166.     int which, i;
  167.     float gr;
  168.  
  169.     which = 0;
  170.     gr = e[which];
  171.     for (i = 1 ; i < 4 ; i++)
  172.     {
  173.         if (fabs(e[i]) > fabs(gr))
  174.         {
  175.             gr = e[i];
  176.             which = i;
  177.         }
  178.     }
  179.  
  180.     e[which] = 0.0;
  181.  
  182.     e[which] = fsqrt(1.0 - (e[0]*e[0] + e[1]*e[1] +
  183.         e[2]*e[2] + e[3]*e[3]));
  184.  
  185.     /* Check to see if we need negative square root */
  186.     if (gr < 0.0)
  187.         e[which] = -e[which];
  188. }
  189.  
  190. /*
  191.  * Build a rotation matrix, given Euler paramaters.
  192.  */
  193. void
  194. build_rotmatrix(Matrix m, float *e)
  195. {
  196.     m[0][0] = 1 - 2.0 * (e[1] * e[1] + e[2] * e[2]);
  197.     m[0][1] = 2.0 * (e[0] * e[1] - e[2] * e[3]);
  198.     m[0][2] = 2.0 * (e[2] * e[0] + e[1] * e[3]);
  199.     m[0][3] = 0.0;
  200.  
  201.     m[1][0] = 2.0 * (e[0] * e[1] + e[2] * e[3]);
  202.     m[1][1] = 1 - 2.0 * (e[2] * e[2] + e[0] * e[0]);
  203.     m[1][2] = 2.0 * (e[1] * e[2] - e[0] * e[3]);
  204.     m[1][3] = 0.0;
  205.  
  206.     m[2][0] = 2.0 * (e[2] * e[0] - e[1] * e[3]);
  207.     m[2][1] = 2.0 * (e[1] * e[2] + e[0] * e[3]);
  208.     m[2][2] = 1 - 2.0 * (e[1] * e[1] + e[0] * e[0]);
  209.     m[2][3] = 0.0;
  210.  
  211.     m[3][0] = 0.0;
  212.     m[3][1] = 0.0;
  213.     m[3][2] = 0.0;
  214.     m[3][3] = 1.0;
  215. }
  216.